% 这是一个例程，研究了多相码(Zadoff-Chu 码 和 Frank 码 )
% 在DAC前后，滤波器以及DAC采样率对包络的影响

clc; clear; close all;

%% 基本参数
N = 64;                   % 码长
oversample_rates = [4, 8, 16, 512];  % 不同过采样率（模拟 DAC 采样率）
fc = 1e6;                 % 基带最大频率
fs_base = 2*fc;         % 最低 Nyquist 采样率
filter_order = 5;

%% 生成 Zadoff-Chu 码（ZC） 和 Frank 码
% Zadoff-Chu 序列 (r relatively prime to N)
r = 5;
n = 0:N-1;
zc_code = exp(-1j*pi*r*n.*(n+1)/N);

% Frank code
L = sqrt(N);
if floor(L) ~= L
    error('Frank码要求码长为完全平方数');
end
frank_code = zeros(1,N);
for m = 0:L-1
    for k = 0:L-1
        frank_code(m*L + k + 1) = exp(1j*2*pi*m*k/L);
    end
end

%% 可视化函数
% plot_envelope = @(x, fs, title_str) ...
%     subplot(2,1,1); plot((0:length(x)-1)/fs, abs(x)); ...
%     title(['Envelope: ' title_str]); xlabel('Time (s)'); ylabel('|x(t)|');

%% 开始处理每个码和每个采样率
codes = {zc_code, frank_code};
code_names = {'Zadoff-Chu', 'Frank'};

for ci = 1:2
    code = codes{ci};
    figure('Name', [code_names{ci} ' Code']);
    
    for si = 1:length(oversample_rates)
        osr = oversample_rates(si);
        fs = osr * fs_base;
        
        % 插值（模拟过采样 DAC 输出）
        code_upsampled = upsample(code, osr);
        code_filtered = conv(code_upsampled, ones(1, osr)/osr, 'same');  % ZOH + smoothing
        
        % 模拟模拟滤波器（巴特沃斯低通）
%         [b, a] = butter(filter_order, 1.2*fc / fs*2);  % 截止频率为fc
        [b, a] = butter(filter_order, 0.5);  % 截止频率为0.25仿真采样率
        code_filtered = filter(b, a, code_filtered);
        
        % 绘制包络
        subplot(length(oversample_rates), 1, si);
        t = (0:length(code_filtered)-1)/fs;
        plot(t*1e6, abs(code_filtered));
        title([code_names{ci} ', OSR = ' num2str(osr) ', Filtered Envelope']);
        xlabel('Time (μs)');
        ylabel('|x(t)|');
    end
end


% 非平滑码
load('non_smoth_weights.mat', 'non_smoth_weights');
len_weight = length(non_smoth_weights);
nonsmoth_code = ifft(non_smoth_weights, 8*len_weight);
% 平滑码
load('smoth_weights_pca_l1.mat', 'smoth_weights');
len_weight = length(smoth_weights);
smoth_code = ifft(smoth_weights, 8*len_weight);


% 随机OFDM
rand_code = ifft(exp(1i * 2 * pi * rand(len_weight, 1)), 4*len_weight);
% 产生LFM采样序列
l = (1:length(smoth_code))';
l = l.^2;
Lfm_sample = exp(1i * pi * l / length(smoth_code) / 8);
% 单频信号，归一化频率为0.2pi
l = (1:length(smoth_code))';
single_tone_sample = exp(1i * 0.15 * pi * l);


codes = {smoth_code,nonsmoth_code,rand_code, Lfm_sample,single_tone_sample};
code_names = {'Smoth', 'Non Smoth','Rand code','Lfm Sample','Single Tone'};

for ci = 1:5
    code = codes{ci};
    figure('Name', [code_names{ci} ' Code']);
    
    for si = 1:length(oversample_rates)
        osr = oversample_rates(si);
        fs = osr * fs_base;
        
        % 插值（模拟过采样 DAC 输出）
        code_upsampled = upsample(code, osr);
        code_filtered = conv(code_upsampled, ones(1, osr), 'same');  % ZOH + smoothing
        
        % 模拟模拟滤波器（巴特沃斯低通）
        [b, a] = butter(filter_order, 1.2*fc / fs*2);  % 截止频率为fc
%         [b, a] = butter(filter_order, 0.25);  % 截止频率为0.25仿真采样率
        code_filtered = filter(b, a, code_filtered);
        
        % 绘制包络
        subplot(length(oversample_rates), 1, si);hold on;legend();
        t = (0:length(code_filtered)-1)/fs;
        plot(t*1e6, abs(code_filtered), 'DisplayName', 'LPF out');
        stem(t*1e6, abs(code_upsampled), 'filled', 'DisplayName', 'upsampled');
        title([code_names{ci} ', OSR = ' num2str(osr) ', Filtered Envelope']);
        xlabel('Time (μs)');
        ylabel('|x(t)|');
    end
end


